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Our study is based on the work of Stinchcombe [J. Phys. C 7, 179 (1974)] and is devoted to the 
calculations of average conductivity of random resistor networks placed on a binary Bethe lattice. 
The binary structure of the Bethe lattice is assumed to represent the normal directions of the regular 
lattice. We calculate the anisotropic conductivity as an expansion in powers of inverse coordination 
number of the Bethe lattice. The expansion terms remained deliver an accurate approximation of 
the conductivity at resistor concentrations above the percolation threshold. This was confirmed by 
means of a Monte Carlo simulation on the regular 2D and 3D lattices. 



I. INTRODUCTION 

The random percolation theory due to Broadbent and 
Hammersley [1] is too simple to explain the great variety 
of percolation phenomena. One confronts complexity of 
real systems with both correlations and anisotropy play- 
ing important role. The motivation for our study is to un- 
derstand better the nature of the anisotropy in electrical 
conductivity of percolating systems. This is approached 
by means of the random resistor network model (RRN) 
[2] . To our knowledge the first analytic solution for con- 
ductivity of anisotropic RRN on the regular square (2D) 
and cubic (3D) lattice was obtained by Bernaskoni [3]. 
Another useful result is that the resistor networks can 
be associated with the networks of saddle points in the 
conductivity profile of high-contrast systems as proved 
in [4]. Currently, anisotropic RRN models are studied 
in geophysics, which could possibly help in prediction of 
earthquakes [5, 6]. Besides conductivity, RRN has been 
used to predict magnetic properties of materials [7] and 
even to estimate sample destruction under critical me- 
chanical stress [8] . The latter research provides a possible 
link to the science of earthquakes too [9] . 

Present study is based on RRN model and uses the 
exact solution on the isotropic Bethe lattice obtained 
by Stinchcombe [10]. This solution was demonstrated 
to describe the macroscopic conductivity curve for RRN 
on the regular 3D lattice calculated with a Monte Carlo 
(MC) simulation [11]. In addition, this solution can be 
used in connection with several important problems in 
solid-state physics [12]. We provide a generalization of 
the Stinchcombe's calculation to the case of a binary 
Bethe lattice, see Fig 1. Besides absence of closed loops, 
this structure has a special feature of being anisotropic at 
each node. Specifically, there are Ua bonds of a kind and 
np bonds of /3 kind connected at each branching point. 
At the same time, their total sum at a node is equal 
to a constant number z refereed to as the coordination 
number of the lattice. We would like to make a distinc- 
tion between the finite Bethe lattice known also as the 
Caley tree and the infinite lattice with the surface sites 
neglected. The latter case is considered in the present 
contribution. 




FIG. 1: Anisotropic Bethe lattice of coordination number z — 
3 with two kinds of bonds, a and /3, depicted by solid and 
dashed lines. The center O, referred to as origin, is where ric 
and np branches are connected by their root bonds of a and 
l3 kind, respectively, z — Ha + np. Each branch is made of 
2 — 1 subbranches connected together by their root bonds. 



The outline of the paper is as follows: in section II we 
present the model and mathematical formulation of the 
problem, in section III we present the main results, in 
section IV we describe our verification of the theory to 
MC simulations, in section V we provide a discussion and 
conclusions. Appendices A-C provide the details of our 
computation. 



II. MODEL 

Unlike several previous anisotropic percolation theo- 
ries [13-17] based on the model of different probabilities 
of filling for two types of lattice bonds, we consider the 
distribution of resistors isotropic. However we make the 
local conductivities of the network elements on the binary 
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Bcthc lattice the 'direction'-dependent, i.e. equal to ctq, 
or Up for a ov (3 occupied bonds, respectively. Thus, the 
resistors are the occupied bonds of the lattice. The lat- 
tice itself is considered non-conductive. Mathematically, 
the local conductivity distribution function is written as 
follows: 



Qaicr) = p5{(T - <Ja) + (1 -p)S{a), 



(1) 



where p is the bond occupation probability, common for 
the bonds of both types. Given p and the conductivities 
of network elements, CTq, and ap, we compute the average 
conductivity of the network connected to a constant po- 
tential source at the origin and grounded at infinity. The 
question how to perform configurational averages turns 
out a difficult one. 

A starting point of the present development is the ob- 
servation that the percolation threshold is given by the 
usual equation [18]: 



Pe= 1/(^-1), 



(2) 



as its derivation does not require considerations of con- 
ductivity as such. This is the consequence of the occupa- 
tion probability common for the bonds of the both kinds. 
The isotropic percolation threshold of sticks of large 
length is expected from the Monte Carlo simulations 
[19] and observed experimentally in polymer/carbon- 
nanotube (CNT) composites [20]. 

Next one defines the probability distribution functions, 
^a{b) and 0/3(6), of the branch conductivities b, such that 



(j)a{b)db= 1. 



(3) 



Those functions measure the contribution from averaging 
over ensemble sampled by resistor permutations, so that 
the average branch conductivity is given by 



f 
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b(t>a{b)db 



(4) 



Here, the symmetry a ^ (3 holds for all quantities. 
Note that we specified in (3) and (4) only the components 
a for the sake of brevity. The second equation is obtained 
readily using a /? interchange. This convention is 
followed everywhere in the text. 

Therefore, the average conductivity of a part of the 
tree consisting of branches connected at the origin in 
parallel is 



ria bo 



(5) 



As an example, one can take ria = n, and np = z — n, 
which leads to n — 1 and z — n oi a- and /3-subbranches, 
respectively, for the a branch (see Fig. 1). In order to 
compute (t)a{b) and (/>/3(6) we use the algorithm of Ref. 
[10] modified to account for the lattice binary structure, 
detailed calculation is placed to the Appendices A-C. 



III. RESULTS 

The analytical solutions have been obtained for the 
two cases: (I) for the case of infinitely large coordination 
number, z — > cx) and (II) near the percolation threshold, 
p ^ Pc- In both cases the solution is represented in the 
form of a Taylor expansion in terms of the small param- 
eters, Pc = {z — 1)^^ and e = {p — Pc)/Pc, respectively. 

In the first case we obtain (see Appendix B for details) 



3c«(I) 



(6) 



where 



= ^As, 



p^ 



G(3) = pA2s[p(2p-l) + 3As], 



= ^Ah[3s^A-2sp 
+A{1 - 3p + 3p^) 
-lOA^/2^-1 



15A3s2 



1 

p2 



= 0(p^), 
A = p — Pc, and s = 1 — p. 



(7) 



This equation gives the average conductivity of a branch 
starting from the a-bond connected to a potential differ- 
ence between the node at its root and the nodes at in- 
finity. The first term is the conductivity of the infinitely 
branched Bethe lattice, while the summation over G-s 
represents the corrections up to and including (z — 1)~^ 
order. Substituting aa = crp into eq. (6) leads to the 
formula that differs from that in Ref. [10] by the fac- 
tor A/p before the summation. We carefully verified this 
and will provide the comparison of our modified and the 
original result to the Monte Carlo (MC) data available 
in literature, see the next Section. 

We now calculate the critical exponents and the 
anisotropy nc;ar the percolation threshold, p ~ Pc- The 
details of the calculation are shown in the Appendix C. 
Close the critical point the integer numbers n/j and 
are set by 



(8) 



received from the symmetry considerations, eq. (C.14). 
Qualitatively, this can be explained as follows: Bethe 
lattice, with its origin O representing a point inside the 
sample, has the branching topology of the infinite cluster. 
Suppose the system is just above pc- The infinite cluster 
should form first in the direction where the resistance to 
current is minimal. For the case when the occupation 
probability p is the same in all directions, the direction 
dependent percolation probability can only be achieved 
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if the fraction of bonds of one kind is larger than another. 
Indeed, the isotropic percolation probability P = 1 — i?^, 
where i? < 1 is the probability to have the finite cluster 
[18], can be generalized to the anisotropic one. Pa — 
1 — which gives Pa > Pp if Ua > np. This leads to 
the following conditions: 



(9) 



Thus, na and np are fixed by the local conductivities. 

Returning for a moment to the previous case, we 
note that the Bethe lattice topology should be intact on 
change of p. Thus, the condition (9) has to be applied 
above the critical point to obtain from cqs. (5) and 
(6): 

In the critical region, we investigate the anisotropy ra- 
tio of the network conductivities and relate this to the 
experimental quantity ct^/oT, where are the bulk 

conductivities parallel and normal to the direction of an 
applied voltage. According to Skal and Shklovskii [22], 

a[\/aT-l^{p-Pcf^'\ (11) 

where \{d) is a critical exponent determined by d - the di- 
mensionality of a problem. Sarychev and Vinogradov [23] 
using the renormalization group theory and computer 
simulations found that A(2) = 0.9 ± 0.1 and A(3) = 
0.3 ±0.1 for 2D and 3D, respectively. Carmona and 
Amarti[24] deduced from experimental data for short car- 
bon fiber reinforced polymers that A(3) ~ 0.4. The de- 
tails of our computation are placed in the Appendix C. 
Our final result (C.21), written in a more concise form, 
is given by 



O'a(II) 



= 0.762- 



^^-^P.,^^0{e% (12) 



2 CTq + (T/3 

where e = {p — Pc)/pc- The discrepancy of the critical 
exponent of 2 with the value of 3 in the case of infinite 
dimensions examined by de Gennes [25], we assign to 
the category of unresolved problems especially because 
de Gennes uses the concept of surface that is never con- 
sidered in Bethe lattice theories. It is known also that 
the critical exponents a and r in the classical Flory- 
Stockmayer theory of the sol-gel transition formulated on 
Bethe lattice are close to those at 3 dimensions [26]. In 
the conductivity field, the critical exponent 2 is assigned 
as the standard one for 3D systems [18, 27]. 

We also calculated the average anisotropy ratio near 
the percolation threshold (C.23): 



C^a(II)/CT/3(jj) - 1 - 



aa(Tfj 



0(6^ 



(13) 



This result is analogous to (11) for the case when a and 
/3 are associated with the parallel and the perpendicular 
components, respectively. Thus, we receive the critical 
exponent, A = 1, which is consistent with the exponent 
obtained in Ref. [28] by an analogous method. 




FIG. 2: Comparison of the near-mean-field formula of Ref. 
[10] (dashed curves) and CT(i)_(n) given by eqs. (10), (12) with 
o'a ~ O-0 (continuous curves) with the MC simulation on 2D 
lattice [21]. The inset D shows a resistor network on 2D lattice 
and the dual network of circles with four nearest neighbors. 
The insets A-C show the networks of circles with multiple 
nearest neighbors. Bethe lattice coordination number used is 
z = 3,5, 10 for the number neighbors on 2D lattice equal to 
4, 8, 16, respectively. 



IV. COMPARISON WITH MONTE CARLO 
SIMULATIONS 

In the work of Strelniker [21] the MC data on a reg- 
ular 2D lattice, for a case when the number of neigh- 
bors is greater than four, are shown. The insets A-D in 
the Fig. 2 demonstrate the model used by Strelniker to 
reduce the resistor network problem on the square lat- 
tice to the problem of multiple nearest neighbors. The 
dashed curves represent the original formula, the contin- 
uous curves show eqs. (10) and (12) with CTq = ap and 
z — 3,5, 10 for the cases A,B,C respectively. We have 
chosen these values to match the data extracted from 
the graph in Ref. [21]. Since the only available limits 
investigated are z — > oo and p ~ Pc and there is no con- 
nection between them, we just glued these two solutions 
at their intersection point as seen in the blow up (Fig. 
2). It is clear that our correction provides a better fit 
well above the percolation threshold, while the original 
formula seems to be smoother close to pc- 

The analytical results (10) and (12) need to be verified 
by the Monte Carlo computer modeling on the regular 
square (2D) and cubic (3D) lattice. For this purpose, 
we performed the simulation of resistor networks with 
the resistor values different for horizontal and vertical 
lattice directions. In the simulation box, two opposite 
boundaries, see 2D scheme in the insets of Fig. 3, are as- 
sociated with the electrodes that are kept at some fixed 
potentials V = const and 0, whereas the other surfaces 
of the box are taken with the periodic boundary condi- 
tions. We used the square lattice of size 50 x 50 and the 
cubic lattice 18 x 18 x 18. Our computer program solves 
the Kirchhoff equations at each lattice node using the 
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FIG. 3: Comparison between the MC simulation on 2D, 3D 
square lattices and the Bethe lattice theory. The theory pre- 
dictions with 2 = 3 and 4 for 2D and 3D, respectively, are 
shown by the curves. The inset shows a resistor network on 
fully occupied regular lattice {p = 1). Resistors are shown by 
two colors to elucidate their orientation-dependent conductiv- 
ities. Two electrodes are connected to an external source of 
constant voltage. Periodic boundary conditions are set in the 
perpendicular direction. 

iteration-relaxation method described in [2]: 1000 ran- 
dom configurations are generated for each p and the a 
parameter from Ref. [2] to get the average conductiv- 
ity. The advantage of the iteration-relaxation method in 
comparison with the alternative transfer-matrix method 
[29] is in relative simplicity of the former. However, the 
lack of the iteration-relaxation method is in progressive 
computational delay upon entering the threshold prox- 
imity. 

Present analytical results and the simulation data are 
compared in Fig. 3. Continuous curves correspond to 
the theory, filled and open squares (circles) represent MC 
data for the average conductivity of occupied bonds of 
type a {(3) on the regular 2D and 3D lattices, respectively. 
For the comparison with these MC data, the Bethe lattice 
of z = 3 and 4 was used for 2D and 3D, respectively. 
Thus, in view of quite good conformity of the theory for 
z = 3, 4 and MC data, the fact of fast convergence of the 
(z — l)~^ expansion is obvious. Furthermore, we find that 
the factor A/p in (6) is crucial for the anisotropic case, 
especially for the averaged component corresponding the 
preferential conductivity direction. 



V. DISCUSSION AND CONCLUSIONS 

We proposed the model of the resistor network that 
has a property of anisotropy in a sense that the conduc- 
tivities of the resistors differ with respect to the lattice 
bond type. We solve the problem in the framework of 
the anisotropic Bethe lattice approximation. The math- 
ematical problem was formulated in terms of a nonlinear 
integral equation, which was solved asymptotically using 



series expansions in two limiting cases: near the percola- 
tion threshold and near the mean- field limit of z —^ oo. 

As noticed in [5], significant anisotropies observed in 
geophysics at the macroscale could be explained by the 
formation of fractal structures in a microscale. Here, the 
macroscopic observable anisotropy is the property of en- 
tire network and the local (intrinsic) anisotropy is asso- 
ciated with the anisotropy in conductivity at a branch- 
ing point of the Bethe lattice. The former is defined 
as a = Oa/aJ}, whereas the latter is essentially the ra- 
tio a ~ Oajafi being the only parameter entering eqs. 
(6,7,10). We find that the present theory is capable of 
producing the strong global anisotropy, a, at small local 
anisotropy, a, in the case when z is large and p 3> Pc- 
Indeed, in this limit the conductivity is given by the 
mean-field formula: o^(i) ~ zp(T^/((TQ-|-cr^), which yields 
0(1) « a} . Interestingly, that previous theories based on 
anisotropic occupation probability [14, 17] predicted the 
opposite: at strong local anisotropy - weak global one. 

Returning to the insets A-C of Fig. 2 we observe that 
the number of the special directions corresponding to the 
vertical or horizontal remains finite as z ^ oo. This is in 
conformity with the assumption z ^ n used in the Ap- 
pendix B, but one has to keep in mind that this is only a 
qualitative explanation. Bethe lattice cannot be embed- 
ded in a finite dimensional space, which suggests that, 
generally speaking, one cannot give a simple geometrical 
picture relating the local anisotropy of the binary Bethe 
lattice to some normal space coordinates. 
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APPENDIX A: ANISOTROPIC BETHE LATTICE 
THEORY [10, 30] 

For a branch and its z — 1 next generation subbranches 
(see Fig. 1) a set of conductivities is defined: 

which are zero or finite according as the corresponding 
root bonds are empty or occupied, where the index i = 
is reserved for the branch origin. These subbranches are 
connected in parallel, so that the conductivities })a^ are 
given by 
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Both ba^ (^/j'') s-iid (Tq, {(7fj), which arc the branch and 
bond conductivities, respectively, are random variables 
distributed with some probability density functions. In 
order to determine the branch distribution functions 
4'a{b^^'^) defined in (3) and (4) we average over various re- 
sistor configurations on the lattice using the distribution 
functions (t>a{b^^^) and 5a (fr*^*^) defined for the conductiv- 
ities of subbranches and individual bonds, respectively, 
so that 



9oc{o-) =p5((T-(Ta) + (l-p)5((T). 



(A.2) 



Being more specific we determine 0q(&) [note that the 
superscript (0) is suppressed for brevity] by performing 
an asymptotic analysis of 



/ pOO POO \ 

/ POO POO \ 



<5(&-6i°)), 



(A.3) 



Since is actually a scries of delta functions, it is 

convenient to introduce the Laplace transform of 4>a{b), 
generally known as the moment-generating function 



B, 



nOC 

= / e-'"'cPa{b)db. 
Jo 



(A.4) 



In the present study, this quantity is named the branch 
generating function. Equation (A.4) combined together 
with (4) leads to 



6a = -5,(0), 



(A.5) 



which means the average branch conductivity is just the 
negative first derivative of Ba{q) evaluated at g = 0. 
Taking the Laplace transform of eq.(A.3) it can be shown 
that 



i?a(g)=a(g)""-lC/3(g)"^ 



where 



CM) 



pOC pOC 

/ daga{cr) / db(/)a{b)exp 
Jo Jo 



qab 
a + b 



(A.6) 



.(A.7) 



Therefore, on account of (A.6), Ca{q) and C/3{q) can be 
named the subbranch generating functions and, in anal- 
ogy with (A.5), one can define the subbranch average 
conductivity as 



b^^ = -c'M- 



(A.8) 



Since 4>a{b) and ga{(^) are the probability densities nor- 
malized to unity, sec (3) and (A.2), respectively, the 
boundary condition for Ca (q) at g = is 



Ca{0) = I. 



(A.9) 



The other boundary condition at g = cxd is identified as 
the probability to have the finite cluster, R, since the 
main contribution to the integral (A.7) comes from the 
neighborhood of 6 = or ct = 0: 



Ca(oo) = R. 



(A.10) 



After some algebra [10] which involves an additional 

Laplace transform that introduces a new variable t, one 
obtains the integral equation 



^ e-"iCc,iq)dq = J dagc,ia){t + a)-^: 



exp 



qat 



Ba{q)dq 



(A.11) 



t + o-Jo \ + 

which is the final exact result to be solved asymptotically. 



APPENDIX B: NEAR-MEAN-FIELD 
EXPANSION 

Consider the integrals on both sides of eq. (A. 11). 
These integrals will be approximated for large t values 
using the Laplace method. [31] The method is based on 
the idea that the main contribution to the integrals comes 
from the neighborhood of q = 0, which makes it possible 
to use the Taylor series expansion as follows 



Ca{q) = e 



ln[Cc(0)-l-gC^(0)-l-...] 



for g<l. (B.l) 



This defines c/a multiplies. Additionally, we have ba^- 
coefficients by 



CaiqY 



„mgC„(0) 



(B.2) 



Substituting (B.2) into (A. 11), one obtains 



f 
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e '''^Ca{q)dq= / duga{u) 



1 



t-Ta{u) {t-\-u)' 



(B.3) 



A;=0 



where 



ut 



Ta{u) = 



--(na-l)CjO)-n^C^(O) (B.4) 
(na-l)C^(O) + n0C:(O) 

(B.5) 



u-(na-l)C^(0)-n^C.(0) 



The formula (B.3) represents an expansion in inverse 
powers of z which is seen from (B.4, B.5). Inversion of 
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the Laplace transform in (B.3) yields 

Ca{q) = j duQaiu) ex.p[qTa{u)] X 

1 



OO k—1 



k=2 



^ (k- r)! ■ 

r=0 ^ ' 



y2(k-r)qk-r 



n-K-l)C^(0)-n;3C7:(0) 



2k- 



where fc-iC'r are the binomial coefficients. The equation 
right above is combined with (B.l), then one equates 
term by term the factors of the successive powers of q, 
and obtains 



= / duQaiu) 



for 



and 



J dug4u)[Taiu) - 

OO 

OO 

E ^i'"^ X 

-71=2 

,fc-l} 

E k>2, 



(B.6) 



m=2 
min{m— 1, k — 1} 



(B.7) 



s=0 



for Z = 1, where 

aWo= r dug^{u) 
Jo 



fc! 



(B.8) 



and 



j{mks) 



-iCs 



(s + l)![fc-(5 + l)]! 



duga(u) X 



2(s+l) 



[t„(z.)-C:(0)]'=-(«+i) 



[M-K-l)C(0)-n/3C;(0)] 



(B.9) 



The first term on the right-hand-side of (B.6) is the rep- 
resentation of the mean- field limit 2; ^ oo, 







duga{u)[Ta{u) - C„(0)], 



(B.IO) 



and the sum over m gives the corrections in inverse pow- 
ers of z. Two equations, obtained by the interchange 
of a and /3 in (B.IO), will be solved neglecting = n 
as it is a constant negligibly small compared to 2; — 1. 
[Here, in order to keep up with the Stinchcome's results, 
wo expand in powers of inverse z — 1 and not z, which is 
equivalent]. Solving (B.IO) we get, as the first solution, 
the isotropic mean-field conductivity: 



Ca(0)i: 



-(Ja{v-Pc)- 



(B.ll) 



Additionally, we obtain 

c'm = - 



(Taf^pip- Pc}V 
(7aPc + 0-^(p-Pc)' 



(B.12) 



which is the anisotropic solution of main interest for us. 

We now move to some elaboration regarding the or- 
ders of correction contained in (B.6)-(B.9). One finds 



that la^^^ba' is of the order (z — 1) ^, since la 

are of the orders {z - l)-('"+i) and {z - l)'"/^, re- 
spectively. Note that the correction to ^° given by 
the first term of the sum in (B.7) affects la^^^ba^ by 
{z — 1)~^ order. Thus, to have the final result up to and 
including 0{[z — 1]""'), the first term in the sum given 

by (B.7) should be taken into account, but only a^^^ 
can be used for to > 2. In addition, when approximat- 
ing aL™' and la'^^^^ with to > 2, we use the replacement 

C^(0) = Cq,(0). This is perfectly acceptable if to > 2, 
since any correction to this would be of {z — 1)~^ or- 



(mlO) 



and 



(m)0 



as (z — 1) 



-2m 



der, and hence would contribute to a, 
To compute the error introduced by this substitution for 

TO = 2, we expand Ta = la^^^aa^^ near Ta = I^^^^^^a^^' 

Ta = f^ + {D"Ta)ACa + {D^'Ta)ACfj, (B.13) 

where AC„ = C^(G) - CJG), D^T^ = 
[5T„/5C^(0)]^/^g^. Additionally, from (B.6) with 

the term m = 2 only, one has 

-fa = {D''Aa)ACa + {D''Ap)AC^ (B.14) 

where D^A^ - [8 J dugaiu)[ra - C:(G)]/5C^(0)]^.^^^. 
The computation of the coefficients yields 

D"Ta = -l + 0([z-l]-i), 

D^Ta = 0, 

D^Ac, = 0, 

D^Aa = {z-l)lT%T- (B.15) 
Equations (B.13)-(B.15) combined together give 

T„ = T„(l + (.- 1)7^6^1 

= f„ + 0([z-l]-^). (B.16) 

Substituting this into (B.6), one obtains the expression 
that contains all corrections up to (2 — 1)~^ order: 







{u) k(«)-C:(0)l+(z-l)ai2)°/r°^x 



r(210) 



1 + (Z- 1)7^0) 1)2^(2)0^(310)^ 



(210) 



I 7-(310) l(3)0 , r(410) 7^(4)0 , r(510) l(5)0 

+7(610) (B.17) 
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Direct computation of the eqs. (B.8) and (B.9) yields 



Jm)0 



-(-ircjor(i~p)x 



i + (-ir 



i-p 



=(mlO) 



j(220) 



C„(0) +(T„p 



m+l 



2(l-p)C„(0) C„(0)+a„p 



(B.18) 

(B.19) 

(B.20) 



Finally, substituting (B.18)- (B.20) into (B.17) and using 
(A. 5), we get the anisotropic conductivity in the form 
of a series expansion in successive powers of the inverse 
coordination number: 



"a (I) = 

where 

G(4) 



(B.21) 



^ A^s [Ss^A - 2sp + A(l - 3p + 3^^) 
+10A2s('^^Vi5AV1 



p ) p- 
p — Pc, and s = 1 — p. 



(B.22) 



The expansion (B.22) coincides with the result of Ref. 
[10] except that factor 5 in front of A^s (^) for G^^^ 
needs to be replaced by 10 according to us. 



APPENDIX C: INVESTIGATION OF CRITICAL 
INDICES 

In this appendix wo investigate the critical exponents 
by means of an asymptotic analysis of the integral equa- 
tion (A. 11) for p approaching pc from above introducing 
a small parameter 



P-Pc 
Pc 



(C.l) 



On the one hand, it is known [18] that the percolation 
probability P, defined on the Bethe lattice as P = 1 — i?^, 
can be expanded near the percolation threshold in series: 

P(e) =Se + Ce2 + 0(e3), 



(C.2) 



where B and C are constants and hence 

ij=l_j(i)e-5(2)e2 + 0(e3). 
Here, 

(5(1) = 2/(z - 2), (C.3) 

while the numerical value of 6^^^ has no significance for 
us, as shown in the analysis set forth below. 

On the other hand, the anisotropic conductivity ex- 
pansion in terms of e has not previously been addressed. 
Motivated by the analysis of Ref. [10], we propose a trial 
solution to (A. 11) of the form 



(C.4) 



(1 2) 

where Ca ' {q) are slowly varying functions of q de- 
scribed by the scaling relations 

Ci'^'Hl) - fL'-'\c^eq), (C.5) 

with Ca being a constant to be determined later. 

Now, the variables s = t/{cae) and y = Casq are de- 
fined. Substituting (C.2) and (C.4) into the left-hand- 
side of (A. 11), multiplying both sides by t and expressing 
the integrals in terms of the new variables, one finds 

(Ta(l + e) , 



/ 

Jo 



where Ca{y) are given by (C.4) and (C.5). A set of equa- 
tions is obtained equating the e-expansion coefficients of 
the same order in the left- and right-hand-side of the in- 
tegral equation (C.6). Firstly, equating the terms linear 
in e we obtain 

fi'Hy)^Pc[{nc, - l)fi'Hy) + np/;\y)]. (C.7) 
As a result, the first correction is isotropic. 

Secondly, equating the terms proportional to and using 
(C.8) yields 



/ e-^yfi^\y)dy = / dye 
Jo Jo 



l + ^{-2s + s^y) 



+ Pc[{no 



(C.9) 

Isotropic solution: When substituting Ca = ctq, we 



recover the isotropic solution /< 



(2) 



f(2) 
J0 



/(^' satisfying 



e-^yf^^\y)dy 



dye" 



{^if^'^-S^'^r + {f^'^-S^'^)[l-2s + s^] 
+Pc[{na-l)f^'Hy)+n0f^'\y)]}. 
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and 



-f 

Jo 



dye' 



-sy 



This gives a simple solution: 



(C.IO) 



(C.ll) 



where ^ is determined by solving numerically the differ- 
ential equation of the second order, 



0, 



(C.12) 



which gives C'(0) = -0.762. 

Anisotropic solution: Here, we try Cq 7^ CTq, in (C.9) to 
obtain another solution. To simplify (C.9), we combine 
it with (C.IO), which yields 



f 

Jo 



e-'yfi^\y)dy = 



(-2s + s2y) 



^"e-«^|(/(i)-5«) 
+Pc[(na - l)fi^\y) + n0f^/\y)]}dy. 
Prom this equation follows that 

(C.13) 

where x is an arbitrary variable. The quantities fa^ (x) 

and (x) have to be symmetric with respect to a ^ /? 
interchange. This condition is satisfied only in the two 
following cases: 



or 



2aaO-/3/{aa + op)- 



(C.14) 
(C.15) 



(C.16) 
(C.17) 



The second set of conditions, eqs. (C.16) and (C.17), is 
unappropriate here due to the condition z » n utilized 
in the Appendix B. Hence, from (C.14) we have 

na= ZGalicFa+CFli), U/s = Z - Ha- (C.18) 

Substitution of (C.18) into (C.13) yields 

fi'H^) - ff\^) = -5(i)a;^"(x)i^^^. (C.19) 

The latter is symmetric with respect to a and [3 inter- 
change. Finally, combining (C. 4, C. 5, C. 3, C.ll, C.15), we 
write the anisotropic solution: 



z-2 



aa(ieq)e + fi^\apeq)e\ (C.20) 



Then, using (A. 5), (5) and (C.18), we obtain the average 
conductivity up to and including order: 

^(11) = -^|0.762-^a,e2-a,/(2)'(g).-^ 

+Ce^ (C.21) 

where C is a constant which could be determined by 
equating the terms proportional to in the expansion 
of eq. (C.6). As a result the first term in the expansion 
(C.21), proportional to e^, is the symmetric one. In ad- 
dition, we can calculate the difference in the derivatives 
of fa^ and f^^ using eqs. (C.19, C.12), which yields 



''a (II) - 



0.762 



^t-Jl{a„-afs)e'. (C.22) 
z — i 



Therefore, we have 



- g;3(ii) _ Z-\ol,-0% 



(C.23) 
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